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Abstract 

We introduce an orbital dependent electron tunneling model and implement it within the atom 
superposition approach for simulating scanning tunneling microscopy (STM) and spectroscopy 
(STS). Applying our computationally efficient method, we analyze the convergence and the orbital 
contributions to the tunneling current and the corrugation of constant current STM images above 
the W(110) surface. In accordance with a previous study, we find corrugation reversals depending 
on the bias voltage and the tip-sample distance. Explaining this effect, we highlight the role of 
the real space shape of the orbitals involved in the tunneling. Moreover, we calculate corrugation 
inversion maps considering different tip models, and find two qualitatively different behaviors based 
on the tip orbital symmetry. From these maps, we predict corrugation reversals for specific tips 
even at positive bias voltages at enlarged tip-sample distances compared to those at negative bias. 

PACS numbers: 68.37.Ef, 71.15.-m, 73.63.-b 
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I. INTRODUCTION 



The experimental use of scanning tunneling microscopy (STM) and spectroscopy (STS) 
has recently gained a boost, particularly for studying magnetic systems [Jfl. The exper- 
imentally observed effects are not staightforward to explain without a proper theoretical 
support {4,5]. One direction of recent theoretical research is focused on extracting surface 
local electronic properties from experimental STS data 6KH|, which is the convolution of 
tip and sample electronic structures. Another research direction is concerned with the simu- 
lation of STM and STS by using different models, mostly based on electronic structure data 



obtained from first principles. Some recently developed approaches highlight t 



of the inclusion of the tip electronic structure into the tunneling model 111 



also possible to separate tip and sample contributions to STS 



The atom superposition approach for simulating STM [14J and STS [15| is a compu 



re importance 
Q, and it is 



13|- 



tationally efficient method. In the present work we go beyond the so-called independent 
orbital approximation and consider a simple model for orbital dependent tunneling. The 
reliability of this new method is demonstrated by the analysis of the tip-sample distance and 
bias voltage dependent corrugation reversal effect observed on the W(110) surface, where 
we find excellent agreement with a previous work 16] at much lower computational costs. 
The computational efficiency of our method enables us to study this effect in more detail. 
We particularly focus on tip effects, and consider ideal tip models with different orbital sym- 
metries, and a more realistic W(110) tip. We find two qualitatively different corrugation 
inversion behaviors based on the tip orbital symmetry. Moreover, we predict corrugation 
reversals for specific tips even at positive bias voltages at enlarged tip-sample distances 
compared to those at negative bias. 

The paper is organized as follows: The theoretical model of the orbital dependent tunnel- 
ing within the atom superposition approach is presented in section [Til Applying this method, 
we investigate the convergence and the orbital contributions of the tunneling current, as well 
as the corrugation reversal of the W(110) surface depending on the applied bias voltage and 
the tip-sample distance in section IIII1 Summary of our findings is found in section IIVI 
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II. THEORETICAL MODEL OF ORBITAL DEPENDENT TUNNELING 



Recently, Palotas et al. developed a three-dimensional (3D) atom superposition approach 



for simulating spin-polarized STM 



(SP-STM) Q and spin-polarized STS (SP-STS) [15| 



based on previous theories The main advance was the inclusion of the tip elec- 

tronic structure, thus, tip dependent tunneling properties can be studied. In the model, it 
is assumed that electrons tunnel through one tip apex atom, and contributions from indi- 
vidual transitions between this apex atom and each of the surface atoms are summed up 
assuming the one-dimensional (ID) Wentzel-Kramers-Brillouin (WKB) approximation for 
electron tunneling processes in all these transitions. The key input is the projected electron 
density of states (PDOS) of the tip apex atom and of the sample surface atoms obtained 
from ab initio electronic structure calculations. In the present paper, first, we review the 
tunneling current and the differential conductance expressions based on the independent or- 
bital approximation for the vacuum decay of electron states, and then we extend this model 
to include a simple orbital dependent tunneling transmission. We consider the non-spin- 
polarized part of the tunneling only, however, this theory can be applied to SP-STM and 
SP-STS in the future. 

Assuming T = K temperature, the tunneling current at the tip position R TIP (x,y, z) 
and at bias voltage V is given by 

r v di 

I(x,y,z,V) = J — (x,y,z,U,V)dU. (1) 
The integrand is the so-called virtual differential conductance, 

£(x,„,,,I7,V)=4 < 2 > 
x T ( E f + ef/ ' ^ V, z))ME T F + eU- eV)n a s (E s F + eU). 

a 

Here, e is the elementary charge, h is the Planck constant, and Ep and E§, are the Fermi 
energies of the tip and the sample surface, respectively. e 2 e 2 /h ensures that the dl/dll is 
correctly measured in the units of A/V. e has been chosen to 1 eV, but its actual value has 
to be determined comparing simulation results to experiments. The sum over a corresponds 
to the atomic superposition and has to be carried out, in principle, over all surface atoms. 
Convergence tests, however, showed that including a relatively small number of atoms in 



the sum provides converged dl/dU values [l2|. T(Ep + eU,V,d a (x,y,z)) is the energy 



and bias dependent tunneling transmission function, which also depends on the distance 
d a (x,y,z) = \Rtip{ x iVi z ) ~ Ra\ between the tip apex and the surface atom labeled by a 
with position vector R^. The tip and sample electronic structures are included into this 
model via projected DOS (PDOS) onto the atoms, i.e., riT{E) and n a s (E) denote projected 
charge DOS onto the tip apex and the ath surface atom, respectively. They can be obtained 
from any suitable electronic structure method. 

Taking the derivative of Eq.([T]) with respect to the bias voltage, the differential conduc- 
tance is obtained. It can be written at the tip position R TIP (x,y, z) and at bias voltage V 
as the sum of three terms, 

^(x, y, z, V) = -^(x, y, z, V, V) + B(x, y, z, V) + D T (x, y, z, V). (3) 

Here, dl/dU has been defined in Eq.flSJ, and B and Dt are the background and tip-derivative 
terms, 

B(x,y,z,V)=e 2 j (4) 
— {E S F + eU, V, d a (x, y, z))n T {E T F + eU - eV)n a s {E s F + eU)dU, 

D T (x,y,z,V) = -e 2 j (5) 
£ T(Ep + eU, V, d a (x, y, z))^jj{E T F + eU - eV)n a s {E s F + eU)dU, 



x 



the trans- 



221, while the 



respectively [15[. The background term, which contains the bias-derivative o: 
mission function, is usually taken into account in recent STS theories ! 
tip- derivative term containing the energy derivative of the tip DOS is rarely considered in 
the recent literature. 

In the spirit of the independent orbital approximation, the transmission probability for 
electrons tunneling between states of atom a on the surface and the tip apex is of the simple 
form, 

T(E S F + eU, V, d a ) = e~ 2K{u ' v)da . (6) 
This corresponds to a spherical exp onential decay of the electron wave functions irrespective 



of their orbital symmetry 20|, |23j, |24j. Assuming the same effective rectangular potential 



barrier between the tip apex and each surface atom, the vacuum decay k can be written as 



K (U,V)=l 2m ( ^±^ -eU), (7) 



where the electron's mass is m, h is the reduced Planck constant, and <ps and 0t are the 
average electron work function of the sample surface and the local electron work function of 
the tip apex, respectively. The method of determining the electron work functions from the 
calculated local electrostatic potential is reported in Ref. 

Next, we extend this tunneling model by taking advantage of the orbital decomposition 
of the electronic structure data and the real space shape of the electron orbitals. The 
PDOS of the sample surface atoms and the tip apex can be decomposed according to orbital 
symmetry, i.e., real spherical harmonics, 7 G {s,p y ,p z ,p x ,d X y,dy Z ,d 3z 2_ r 2,d xz ,d x 2_ y 2}, as 

n%{E) = Y,nU E )i (8) 

/s 

n T (E) = $> T7 (£). (9) 

7 

Assuming such an orbital decomposition, the virtual differential conductance can be gener- 
alized as 

§ {x ,y, z ,U,V)=^ (10) 
x J2J2 T ^ E f + eU > V > f> z )) n Ti(E T F + eU- eV)n a S p{E s F + eU), 

a 0,7 

where, additionally to the atomic superposition (sum over a), we sum up each (3 -H- 7 
transition with an orbital dependent tunneling transmission: T^(E, V, d a ) gives the electron 
tunneling probability at energy E from the tip apex 7 orbital to the orbital of the ath 
surface atom at positive bias voltage (V > 0), and from the (3 orbital of the ath surface 
atom to the tip apex 7 orbital at negative bias (V < 0). can be defined in different ways 
based on physical arguments. We consider the following form, 

T M {E S F + eU, V, d a ) = e- 2 « u W d "tp 7 (# a , <p a ) (11) 

for each surface atom -B- tip apex ID electron transition. Here, the energy and bias dependent 
part corresponds to the spherical exponential decay considered in Eq.(jSJ), and is independent 
of the orbital symmetry. This is multiplied by an orbital dependent expression tp y , which 
depends on the spatial arrangement of the sample atoms relative to the tip apex and all 
the orbital shapes involved in the tunneling. The angular dependence on § a and y? a comes 
into play in the following way: Let us consider one transition between surface atom a and 
the tip apex along the line connecting the two atoms. A particular geometry is shown in 



Figure [TJ For brevity, we omitted the a notation of the surface atom. For both atoms, a 
local coordinate system can be set up, and the angular dependence of the atomic orbital 
wave functions, (p), are defined in the corresponding coordinate system, as summarized 
in Table [H R TIP — R^ defines a vector pointing from the surface atom toward the tip apex, 
and it can be represented by the (d a ,$ a ,(pa) coordinates in the spherical coordinate system 
centered on the ath surface atom. Taking the opposite connecting vector from the tip apex 
toward the surface atom, its coordinates are (d a , ir — i? , ir + <p a ) in the spherical coordinate 
system centered on the tip apex. According to Figure [TJ d a , i? OJ and (p a can be calculated as 



using the global coordinates R TIP = (x,y,z) and R a = (x a ,y a ,z a ). Considering above, tg 7 
accounts for the modification of the perfect spherical exponential decay along the connecting 
line through the angular dependence of the atomic orbitals as 



where xp,-y{$, V 5 ) are the real spherical harmonics summarized in Table [B They were chosen in 
such a way that < tg 7 < 1. This factor takes the effect of the directional tunneling between 
real space orbitals into account. The physical motivation is the angular dependence of the 
electron density of the orbitals in the given tunneling direction, which modifies the tunneling 
transmission. The maximum t — 1 is obtained if the angular distributions of the electron 
density according to the given orbital symmetries ((3, 7) around both the sample surface atom 
and the tip apex have maxima along the line of the tunneling direction. This is always the 
case for s-s type of tunneling irrespective of the relative position of the tip apex and sample 
surface atoms, i.e., we observe perfect spherical exponential decay between tip and sample s 
orbitals. In some particular geometries t — 1 can be obtained even for other type of orbitals, 
e.g., if the tip apex is precisely above surface atom a, i.e., $ a = 0, then t/3 7 (^ a = 0, ip a ) = 1 
for all of the following combinations: 0, 7 G {s,p z , d 3z 2_ r 2}. On the other hand, if the tip 
apex is above surface atom a, then orbitals with nodal planes orthogonal to the surface have 
zero contribution to the tunneling from this particular surface atom, i.e., a reduced effective 




(12) 
(13) 



(14) 



(15) 
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tunneling transmission is obtained [25|. Note that the independent orbital approximation 
corresponds to tg 7 = 1 for all /3 <H- 7 transitions, i.e., the same tunneling transmission 
is assumed between all orbitals. Within our orbital dependent tunneling approach, ideal 
tip models with particular orbital symmetries can be considered, i.e., 70 orbital symmetry 
corresponds to the choice of nr 70 = l(e^) -1 and nr( 7 ^ 7o ) = 0. More realistic tips can be 
obtained by explicitly calculating the orbital decomposition of the tip apex PDOS in model 
tip geometries, e.g., in the present work a model W(110) tip is used. 

Our theory is, thus, an extension of the atom superposition STM/STS approach con- 
sidering tunneling between directional orbitals. The tunneling current and the differential 
conductance can be calculated at the tip position R TIP (x, y, z) and at bias voltage V as the 
sum of all (3 -H- 7 transitions between real space orbitals, 

I(x,y,z, V) = ^2,I M (x,y, z, V), (16) 

/S,7 

§(x,y,z,v) = j2^r( x >y> z >n (17) 



respectively, with 



dy y >*> > j ^ dV 

0,-y 



,2 



I 01 (x,y,z,V) =e 2 - (18) 
V / T M {E S F + eU, V, d a (x, y, z))n Tl {E T F + eU - eV)n a sp {E s F + eU)dU, 

a J 

^(x,„,,,10 = 4 (19) 

£ 7> 7 (#f + eV, V, d a (x, y, z))n Tl {E T F )n a s ^E s F + eV) 

. a 

+ J2 / ^(E% + eU,V,d a (x,y,z))n Tl (ET + eU-eV)n a S p(E§, + eU)dU 

- t ^ e f + eU > v > d °( x > y> Z ))^f( E F + eU ~ eV ) n U E F + eU ) du 

This decomposition enables the analysis of the orbital contributions to the total tunneling 



current and to the differential conductance. In relation to Chen's derivative rule [26|, we 
can state that while it is formulated inspired by the Tersoff-Hamann model, and calculates 
the tunneling transmission as the absolute value square of the tunneling matrix element, 
that is proportional to the sample wave function derivative with respect to the real space 



7 



coordinate corresponding to the given tip orbital symmetry (7), our transmission function 
also depends on the sample surface atoms' orbital symmetry /3. Moreover, in our theory 
the electronic structure of the tip apex via the PDOS is accounted for. Note that the 
presented method can also be applied to magnetic systems taking into account the orbital- 



lol together with the orbital 



decomposed magnetization PDOS of the tip and sample fl^ . 
dependent tunneling transmission in Eq. fllll) . As it was pointed out by Ferriani et al. [27|, 
the spin polarization in the vacuum can have an opposite sign than within the tip apex 
atom, and this sign change is also accompanied by different dominating orbital characters. 
Thus, the consideration of an orbital dependent tunneling transmission might be a better 
model for describing electron transport through a magnetic tunnel junction. We return to 
the related spin-polarized STM/STS model in the future. 

III. RESULTS AND DISCUSSION 

In order to demonstrate the reliability of our orbital dependent tunneling model we 
consider a W(110) surface. This surface has technological importance as it is widely used as 



substrate for thin film growth, see e.g., Refs. |16l. I28j. As it was pointed out by Heinze et al. 
[hI , the determination of the position of surface atomic sites is not straightforward, as atomic 
resolution is lost at negative bias voltages, and a bias-dependent corrugation reversal has 
been predicted. This means that normal and anticorrugated constant current STM images 
can be obtained in certain bias voltage ranges, and the W atoms do not always appear as 



protrusions in the images. It was shown that a competition 
parts of the surface Brillouin zone is responsible for this effect 



between states from different 



16. 



29| | . We reinvestigate this 



corrugation reversal effect, as it provides a challenge for our orbital dependent tunneling 



model. We find excellent agreement with the results of Ref. [16J, where an s-wave tip has 
been used, and since our method is computationally more efficient, we study this effect in 
more detail. We particularly focus on tip effects, and consider ideal tip models with different 
orbital symmetries, and a more realistic W(110) tip. We find two qualitatively different 
corrugation inversion behaviors based on the tip orbital symmetry. Moreover, we predict 
corrugation reversals for specific tips even at positive bias voltages at enlarged tip-sample 
distances compared to those at negative bias. Explaining our findings, we highlight the 
role of real space orbital orientational overlaps between the surface and the tip rather than 
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considering electron states in the reciprocal space, thus, a different kind of understanding is 
provided. 

We performed geometry relaxation and electronic structure calculations based on the den- 



pproximation (GGA) imple- 



sity functional theory (DFT) within the generalized gradient a 
mented in the Vienna Ab-initio Simulation Package (VASP) 30 
for electronic wave function expansion together with the projector augmented wave (PAW) 



32| . A plane wave basis set 



method [33j has been applied, and the exchange-correlation functional is parametrized ac- 
cording to Perdew and Wang (PW91) 34 1. 

We model the W(110) surface by a slab of nine layers, where the two topmost W layers 
have been fully relaxed. After relaxation the W-W interlayer distance between the two 
topmost layers is reduced by 3.3%, while the underneath W-W interlayer distance increased 
by 1.1% compared to bulk W. A separating vacuum region of 18 A width in the surface 
normal (z) direction has been set up between neighboring supercell slabs. The average 
electron work function above the surface is calculated to be 4>s = 4.8 eV. We used al5xl5x3 
Monkhorst-Pack (MP) [sjj] k-point grid for calculating the orbital-decomposed projected 
electron DOS, ngp(E), onto the surface W atom. The unit cell of the W(110) surface 
(shaded area) and the rectangular scan area for the tunneling current simulation are shown 
in Figure [2j In our calculations we used the experimental lattice constant a = 316.52 pm. 
Moreover, the surface top (T) and hollow (H) positions are explicitly shown. 

We considered different tip models. The orbital-independent ideal tip is characterized 
by tfr = 1 and n Tl {E) = l/9(eV) _1 , so that n T {E) = E 7 ri T 7 (^) = l{eV)~ l . This ideal 
electronically flat tip represents the limiting case of the independent orbital approximation 
used in previous atom superposition tunneling models 12, 14, 15, 2o|. In order to study the 



effect of the orbital dependent tunneling, other tip models are needed. First, we consider 
ideal tip models having a particular orbital symmetry 7 . In this case, tg 7 is calculated 
following Eq.( fT5"l) . and for the energy dependence of the tip PDOS, nx l0 = l(e^) -1 and 
rir( 7 ^ 70 ) = are assumed. More realistic tips can also be employed by calculating the 
orbital decomposition of the tip apex PDOS in model tip geometries, and using Eq. (1151) for 
the orbital dependent transmission factor. We used a W(110) tip. It has been modeled by 
a slab consisting of seven atomic layers, having one W apex atom on both surfaces, i.e., 
with a double vacuum boundary. In this system the apex atom and the topmost surface 
layers have been relaxed on both sides. The interlayer distance between the topmost W 
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layer and the underneath one was reduced by 0.9%, while the adatom-topmost layer vertical 
distance decreased by 21.2% compared to bulk W. The interaction between apex atoms in 
neighboring supercells in the lateral direction is minimized by choosing a 3 x 3 surface cell, 
and a 16.7 A wide separating vacuum region in the z direction. Moreover, an llxllxl MP 
k-point grid has been chosen for obtaining the orbital-decomposed projected DOS, tit^E), 
onto the apex atom. The local electron work function above the tip apex is 0t = 4.36 eV. 

Previously, the convergence of the dl/dll part of the differential conductance has been 
investigated with respect to the number of surface atoms involved in the summation of the 



orbital-independent atomic superposition formula [12j. Due to the spherical exponential 
decay assumed for the electron wave functions, a rapid convergence was found. We report a 
similar convergence test for the orbital- dependent tunneling approach, comparing different 
tip models. In order to take into account a wide energy range around the Fermi level, we 
calculated the tunneling current at -2.5 V and +2.5 V bias voltages at z = 4.5 A above 
a surface W atom, and averaged these current values. We considered ideal tips of the 
orbital-independent model, and with s, p z , and d^ z 2_ r 2 symmetry, as well as the W(110) 
tip. In order to obtain comparable results, we normalized the averaged current for each tip 
calculation. The convergences of the normalized averaged current with respect to the lateral 
distance on the surface, d\\, characteristic for the number of atoms involved in the atomic 
superposition, are shown in Figure [3) By calculating the current, contributions from surface 
atoms within a radius of du measured from the W atom below the tip apex are summed up 
(sum over a). We find that the orbital-independent, the s-type, and the W(110) tips behave 
quite similarly concerning the current convergence, while for the p z - and c^^-type tips a 
faster convergence is found. This rapid convergence can be explained by the more localized 
character of the latter tip orbitals in the direction normal to the sample surface (z). On the 
other hand, the orbital-independent tip with T = e~ 2Kd is a good approximation for the s- 
type tip (with index 7 = 1), where the spherically decaying transmission function part is still 
dominant, i.e., Tp t i = e~ 2Kd x 2 p because Xi — 1- I n case °f the W(110) tip, electronic states 
of all considered symmetries have a contribution, and their relative importance is not only 
determined by the transmission function via the orbital shapes, but also by the product of 
the symmetry-decomposed electron PDOS of the surface and the tip. In general, the orbitals 
localized in different than the z direction can show a slower current convergence than the 
s orbitals. However, the partial PDOS of such states is relatively low, and, interestingly, 
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we obtain a similar current convergence in the studied energy range as for the s-type tip. 
Choosing different bias voltages, and, thus, different electon states involved in the tunneling, 
for the W(110) tip we found current convergences dissimilar to the s-type tip behavior. The 
convergence can be slower or faster than obtained for the s-type tip, depending on the partial 
PDOS of each directional orbitals in the given energy range. 

Based on the convergence tests, atom contributions from at least d\\ = 3a m 9.5 A distance 
from the surface-projected tip position shall be considered. In case of calculating STM 
images, d\\ = 3a ~ 9.5 A has to be measured from the edge of the scan area in all directions 
in order to avoid distortion of the image, thus, involving 67 surface atoms in the atomic 
superposition. For brevity, in the following, we use the same surface atoms to calculate 
single point tunneling properties as well. 

Let us analyze the relative importance of all /3 •f-)- 7 transitions in determining the total 
tunneling current at different tip positions. From this analysis we obtain a qualitative 
picture about the role of the different atomic orbitals in the construction of the tunneling 
current. The 7g 7 current contributions can be calculated according to Eq. (118j) . These can be 
represented by a current histogram, that gives the percentual contributions of the individual 
transitions to the total current. Figure H] shows such histograms using the W(110) tip at 
V— -0.1 V bias voltage, z = 4.5 A above two different tip positions: part a) corresponds to 
the tip apex above the surface top position, and part b) to the tip apex above the surface 
hollow position, T and H in Figure [2} respectively. We obtain a 9 x 9 matrix from the 
considered orbitals, which are denoted by numbers 1 to 9 following the indices reported 
in Table HI We find that most contributions are due to the s (1), p z (3), d yz (6), d 3z 2_ T 2 
(7), and d xz (8) orbitals and their combinations. The largest contribution to the current 
is given by the d 3z 2_ r 2 — d 3z 2_ r 2 (7-7) transition in both considered tip positions. While it 
corresponds to 25 per cent above the top position, it is slightly lower, 18 per cent above 
the hollow position. Concomitantly, above the hollow position, the relative importance of 
both tip and sample d yz (6) and d xz (8) orbitals is increased as it is expected from the 
geometrical setup, i.e., the d yz — d yz (6-6), d yz — d 3z 2_ r 2 (6-7 and 7-6), d xz — d xz (8-8), and 
d xz — d 3z 2_ r 2 (7-8 and 8-7) contributions correspond to larger orientational overlap of the 
mentioned tip and sample orbitals if the tip is above the hollow position rather than above 
the top position, as suggested by Fig. [2] and Eq. ffToT) . Thus, our simple orbital dependent 
tunneling model captures the effect of the localized orbitals and goes beyond the spherical 



11 



Tersoff-Hamann model. Note that if a larger bias voltage is considered, i.e., the electronic 
states are somewhat averaged, then the independent orbital approach might turn out to be 



a good approximation [20]. 



The role of the localized orbitals can best be demonstrated by returning to the corrugation 



inversion phenomenon found, e.g., on (110) transition metal surfaces [16]. Under certain 
circumstances the apparent height of W atoms at the surface top position (zt) can be larger 
or smaller than the apparent height of the surface hollow position (zh) at constant current 
(/ = const) condition. (For the surface top (T) and hollow (H) positions, see Figure [21) 
Thus, the sign change of Az(I) = zt{I) — zh{I) is indicative for the corrugation inversion. 
Obviously, Az(I) > corresponds to a normal STM image, where the W atoms appear 
as protrusions, and Az(I) < to an anticorrugated image. Since the tunneling current is 
monotonously decreasing with the increasing tip-sample distance, we can obtain information 
about the occurrence of the corrugation inversion simply by calculating the current difference 
between tip positions above the top and hollow sites of the W(110) surface. The current 
difference at tip-sample distance z and at bias voltage V is defined as 

AI(z, V) = I T (z, V) - I H (z, V). (20) 

This quantity can be calculated for specific tips, and we call the AI(z, V) — contour as 
the corrugation inversion map. This gives the (z, V) combinations where the corrugation 
inversion occurs. The sign of AI(z, V) provides the corrugation character of an STM image 
in the given (z, V) regime. Due to the monotonously decreasing character of the tunneling 
current, AI(z,V) > corresponds to Az(I(V)) > 0, i.e., normal corrugation, and the 
same holds for the opposite relation, i.e., AI(z,V) < corresponds to Az(I(V)) < and 
anticorrugation. 

First, we calculated AI(z,V) using the independent orbital approximation and Eq.(|6]) 
for the tunneling transmission, and found that AI(z, V) is always positive. This means that 
the spherical exponential decay itself can not account for the observed corrugation inversion 
effect, and the W atoms always appear as protrusions in STM images calculated with this 
model. However, considering the orbital dependent tunneling transmission, Eq.l flTT) . we find 
evidence for the corrugation inversion effect, thus, highlighting the role of the real space 
shape of electron orbitals involved in the tunneling. Figure shows AI(z, V) = contours 
calculated with different tip models in the [0 A, 14 A] tip-sample distance and [-2 V,+2 V] 
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ipreviously not 



16] , where an s- 



bias voltage range. Before turning to the analysis of the results obtained wit 
considered tip models, let us compare our results with those of Heinze et al. 
wave tip model has been used. They found corrugation reversal at around -0.4 V at z = 4.6 
A tip-sample distance, and above that voltage normal, while below that anticorrugated 
STM images were obtained. Our model with an s-tip provides the same type of corrugation 
reversal at -0.23 V at the same distance as can be seen in part b) of Figure [5] (curve with 
square symbol). These values are in reasonable agreement, particularly concerning their 
negative bias value, where atomic resolution is difficult to achieve experimentally, which is 



an indication for being close to the corrugation inversion regime [16] . On the other hand, 
a linear dependence of the corrugation reversal voltage and the tip-sample distance has 
been reported by Heinze et al. (z = 4.6 A, V = -0.4 V) to (z = 7.2 k,V = V). Our 
model qualitatively reproduces this linear dependence in the same bias range, though the 
quantitative values are somewhat different. 

Calculating the corrugation inversion maps with more tip models, we find two distinct 
behaviors depending on the tip orbital symmetry. Parts a) and b) of Figure [5] show these. 
While, there is no corrugation inversion occurring at positive bias voltages for tip models in 
part a), tips in part b) show this effect also at positive bias. Moreover, in both parts anti- 
corrugation is observed in the large tip-sample distance reg ion (z > 13.5 A) as AI(z, V) < 



there. This is in accordance with the prediction of Ref. [29] based on the analysis of the com- 
peting electron states in the surface Brillouin zone for an Fe(001) surface. In the z < 13.5 
A range, however, the corrugation character in the two parts of Figure [5] is remarkably dif- 
ferent. The tip models in part a) always show anticorrugation at positive bias voltages, 
and below -0.27 V they provide corrugation characters starting from anticorrugation, then 
normal corrugation, and again anticorrugation with increasing tip-sample distance. On the 
other hand, in part b), normal corrugation is found close to the surface, which reverts only 
once with increasing tip-sample distance in the full studied bias range. These different be- 
haviors can be attributed to the tip orbital characters. It is interesting to notice that none 
of the considered tip orbitals in part a) are localized in the ^-direction, and they have nodal 
planes either in the yz plane (p x and d xz ) or in the xz plane (p y and d yz ) or in the x = y and 
x = —y planes (d x 2_ y 2). On the other hand, in part b) there are tips which are localized in 
the z-direction (p z and d^ z 2_ r 2) or having nodal planes in both the xz and yz planes (d xy ) as 
well as the spherical s tip and the W(110) tip that contains all type of orbitals with energy 



13 



dependent partial PDOS functions. The particular tip nodal planes restrict the collection 
of surface atom contributions to specific regions on the sample surface. Furthermore, by 
changing the tip-sample distance, the orientational overlaps between tip and sample orbitals 
change, and according to our model some localized orbitals gain more importance in the 
tunneling contribution, see also Figure SJ Since we calculate the current difference between 
tip positions above the surface top and hollow sites, the complex tip-sample and bias volt- 
age dependent effect of the real space orbitals on the tunneling can be visualized via the 
corrugation inversion maps. 



Concerning tips with p z and d 3z 2_ r 2 orbital symmetry, Heinze et al. [16| calculated a cor- 
rugation enhancement factor of 2 and 6.25, respectively, based on Chen's derivative rule 26]. 
Moreover, they argued that the corrugation inversion map should be practically identical to 
the one obtained by using the s-tip model, and the corrugation values just have to be scaled 
up by these factors. On the contrary, based on our orbital dependent tunneling model, we 
find that the p z and d% z 2_ r 2 tips provide qualitatively different corrugation inversion maps, 
i.e., although their bias dependent shape is similar to the one of the s-tip, their tip-sample 
distance is systematically pushed to larger values, see part b) of Figure El This is due to 
the more localized character of these tip orbitals in the ^-direction. 

Corrugation inversion with the d xy tip occurs at the largest tip-sample distance. A 
possible explanation can be based on its xz and yz nodal planes. While above the top 
position only the underlying W atom, above the hollow position all four nearest neighbor 
W atoms give zero contribution to the current, thus, It is expected to be higher than Ih 
at small tip-sample distances. To overcome this effect, the tip has to move far from the 
surface, since then the relative importance of the nearest neighbor contributions decays 
rapidly compared to other parts of the surface. 

For a more realistic W(110) tip, we find that the corrugation inversion map is qualitatively 
similar to the s-tip. In the [-0.3 V,-0.15 V] bias range they are even almost identical, and 
the same linear dependence is observed. The similarity of the W(110) tip to the s-tip has 
also been found for the current convergence in Figure [31 and it suggests that the common 
effect of different orbitals of a realistic tip can indeed be close to a spherical tip behavior. 
However, we stress again that the independent orbital approximation could not reproduce 
the corrugation inversion effect. 

Apart from above findings, we obtain corrugation inversion also in the positive bias range 
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at enlarged tip-sample distances for the s, p z , d 3z 2_ r 2, and W(110) tips considered in part 
b) of Figure |5j This is most probably due to the surface electronic structure and this effect 
is even more difficult to capture in experiments as the corrugation values themselves decay 
rapidly with increasing tip-sample distance. 

In order to demonstrate the corrugation inversion more apparently, constant current STM 
images can be simulated. As it is clear from Figure[5l any type of crossing of the AI(z, V) = 
contour results in the occurrence of the corrugation reversal. In experiments, two ways can 
be considered to record STM images in the normal and anticorrugated regimes: 1) keep the 
tip-sample distance z constant, and change the bias voltage V, or 2) keep the bias voltage 
V constant and change the tip-sample distance. These modes correspond respectively to a 
horizontal and a vertical crossing of the AI(z, V) — contour in the (z, V) plane as shown in 
Figure |5j As the second option seems to be more feasible experimentally, we simulated STM 
images at a fixed bias voltage of -0.25 V using the W(110) tip, and extracted three different 
constant current contours. As part b) of Figure |5] shows, at this particular bias voltage a 
corrugation inversion at z = 4 A occurs. Therefore, we chose current values that correspond 
to current contours with three different apparent heights of the W atom at zt = 3.5 A, 
4 A, and 4.5 A. We calculated the tunneling current in a box above the rectangular scan 
area shown in Figure [2] containing 99000 (30 x 22 x 150) grid points with a 0.15 A lateral 
and 0.05 A vertical resolution. The constant current contours are extracted following the 
method described in Ref. [ljj]. The obtained STM images are shown in Figure [6] In part a) 
Az(I) = Zt(I) — Zh{I) > 0, thus, the apparent height of the W atom is larger than the one 
of the hollow position. The reverse situation is obtained in part c). Part b) corresponds to 
the situation where the corrugation inversion occurs, i.e., Az(I) = 0. In this case we obtain 
a stripe- like image. This, and all of the simulated STM images are in good agreement with 



Ref. 



16j . We also calculated the corrugation of the individual current contours, and obtained 



0.8, 0.3, and 0.2 pm for parts a), b), and c), respectively. These values are extremely small, 
but are again in good agreement with the reported values in Ref. 16 ]. 

Thus, employing our new orbital dependent tunneling model, we could reproduce and 
reinvestigate the corrugation inversion effect observed on the W(110) surface. Although 
this effect is driven by the surface electronic structure, we showed that different tips can 
drastically modify its tip-sample distance and bias voltage dependence. 
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IV. CONCLUSIONS 



We developed an orbital dependent electron tunneling model and implemented it within 
the atom superposition approach for simulating STM and STS. Applying our computa- 
tionally efficient method, we analyzed the convergence and the orbital contributions to the 
tunneling current above the W(110) surface. We found that the d^ z 2_ r 2 — dj, z 2_ r 2 contribution 
is the largest, and depending on the tip position, other d states can gain importance as well. 
We also studied the corrugation inversion effect. Using the independent orbital approxi- 
mation, no corrugation reversal has been obtained at all. Employing the orbital dependent 
model, we found corrugation reversals depending on the bias voltage and the tip-sample 
distance, in accordance with the work of Heinze et al. 16j. Explaining this effect, we high- 
lighted the role of the real space shape of the orbitals involved in the tunneling. Moreover, 
we calculated corrugation inversion maps considering different tip models, and found two 
qualitatively different behaviors based on the tip orbital symmetry. From these maps, we 
predicted corrugation reversals for specific tips even at positive bias voltages at enlarged 
tip-sample distances compared to those at negative bias. Simulation of STM images made 
the corrugation inversion effect more apparent. Extending this orbital dependent tunneling 
model to magnetic junctions is expected to enable the study of the interplay of real space 
orbital and spin polarization effects in SP-STM and SP-STS experiments in the future. 
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TABLE I: Real space orbitals, their indices used in the present paper, their definition from spherical 
harmonics YJ m ($, <p), and the angular dependence of their wave functions x(i?, (p). Note that i? and 
p are the usual polar and azimuthal angles, respectively, in the spherical coordinate system centered 
on the corresponding (tip or sample) atom. 



orbital 


index 


definition 




s 


1 


Y° 


1 


Fy 
Pz 
Px 


2 
3 
4 


Y} - Yr 1 
i i 

Y^+rr 1 


sin($)sin(p) 

cosifl) 
sin{$)cos((p) 


dxy 
dyz 
d 3z 2_ r 2 

dxz 


5 
6 
7 
8 
9 


Yi - Kf 2 
Yi ~ Y^ 

Y 2 ° 
Yi + Y,' 1 

Yi + y 2 - 2 


sin 2 ('d)sin(2<p) 
sin(2 , d)sin(p) 

\ (3cos 2 (tf) - 1) 
sin{2-ff)cos{ip) 

sin 2 ('&)cos(2p) 




FIG. 1: Geometry of a general tip-sample setup. The a notation of the surface atom is omitted for 
brevity. 
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FIG. 2: The surface unit cell of W(110) (shaded area) and the rectangular scan area for the 
tunneling current simulations. Circles denote the W atoms. The top (T) and hollow (H) positions 
are explicitly shown. 
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FIG. 3: (Color online) Convergence of the normalized averaged current calculated with different 
tip models. 
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FIG. 4: (Color online) Histograms of the current contributions from all tip-sample transitions 
with different orbital symmetries using the W(110) tip at V= -0.1 V bias, a) Tip apex z = 4.5 
A above the surface top (T) position (W atom); b) tip apex z = 4.5 A above the surface hollow (H) 
position, see also Figure [2j The indices of the atomic orbitals (1-9) follow the notation reported in 
Table E 




FIG. 5: (Color online) The AI(z,V) = It(z,V) — Ih(z,V) = contours indicative for corruga- 
tion inversion calculated with different tip models, a) and b) parts show two distinct behaviors 
depending on the orbital symmetry of the tips. The large tip-sample distance region (z > 13.5 A) 
corresponds to AI < in both a) and b) parts. 
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FIG. 6: (Color online) Simulated STM images (constant current contours) at V= -0.25 V bias 
voltage using the W(110) tip. The scan area corresponds to the rectangle shown in Figure (2) Light 
and dark areas denote higher and lower apparent heights, respectively. The current values (I), the 
apparent heights of the W atom (zt), and the corrugation of the contours (Az') are as follows: a) 
I = 2.5 X 1CT 1 nA, z T = 3.5 A, Az' = 0.8 pm; b) I = 9.3 x 10~ 2 nA, z T = 4.0 A, Az' = 0.3 pm; 
c) / = 3.4 x 10~ 2 nA, z T = 4.5 A, Az' = 0.2 pm. 
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